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FlowPy is a numerical toolbox for the solution of partial differential equations encountered in Func- 
tional Renormalization Group equations. This toolbox compiles flow equations to fast machine 
" ■ code and is able to handle coupled systems of flow equations with full momentum dependence, 

I ■ which furthermore may be given implicitly. 

o" 
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O ■ I. INTRODUCTION 

t , The functional renormalization group (FRG) has in recent years been successfully applied to a wide variety of non- 

■ perturbativeproblems such as critical phenomena, fermionic systems, gauge theories, supersymmetry and quantum 
gravity, see |ll4l2| for reviews. Although these systems are very different in their physical nature, the flow equations 

■ always have a similar structure. 

The aim of the FRG is to calculate the generating functionals of IPI correlation functions from which the dynamics 
of the theory can be inferred. The core ingredient is the scale dependent effective action denoted by Fj, with the RG 
scale k. It interpolates between a microscopic description through the classical action at some UV scale k = A and a 
macroscopic description at low energy scales fc = through the full quantum effective action. The RG scale k serves 
as an infrared regulator suppressing all fluctuations with momentum smaller than k. Thus, for fc = all fluctuations 
are taken into account and we have obtained a full solution of the quantum theory. The flow of the scale dependent 
Q ■ effective action is governed by the Wetterich equation ^] 
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with FJ,, being the second functional derivative of the effective action. The momentum-dependent regulator function 
Rk in Eq. ([ij establishes the IR suppression of modes below fc. In the general case, three properties of the regulator 
Rk{p) are essential: (i) i?fe(p)|p2/j,2_j.g > which implements the IR regularization, (ii) Rk{p)\k'^ /p^^o = which 
guarantees that the regulator vanishes for fc — > 0, (iii) Rk{p)\k^K^oo oo which serves to fix the theory at the 
' classical action in the UV. Different functional forms of Rk correspond to different RG trajectories manifesting the 
^\ • RG scheme dependence, but the end point F^^o — > T remains invariant. 

' Solving this partial nonlinear differential equation (JlJ head on is impossible in most cases. Thus approximations 
. for the effective action have to be introduced resulting in a system of coupled differential equations. Recently, a 
■ Mathematica extension which is able to derive Dyson- Schwinger equations and functional renormalization group 



psj . equations was published [1^. However, solving these systems beyond the most simple approximations is a numerical 
T-H ■ challenge. For some systems, e. g. supersymmetric quantum mechanics [isj and the two dimensional A/" = 1 Wess- 
^ I Zumino model the coefficient function of the highest derivative can become singular. This implies that solving the 
• ^ ■ differential equations numerically has to be done with great care. If the equations are solved exactly, the singularity 
■ is never reached but rather the fiow is repelled if it comes close to the singularity. A numerical solution has to take 
^ . this behavior into account. 
■ - - ■ Including a full momentum dependence has become more important over the last few years. This leads to a much 
higher numerical effort for the solution of the flow equations. Full momentum dependence of propagators and vertices 
has previously been treated successfully in the literature [l6l - [23 |. 

This article presents a numerical toolbox called FlowPy for the solution of a broad class of partial differential 
equations that are encountered in the study of the Functional Renormalization Group equations. Specifically, FlowPy 
supports full momentum dependence of the flowing function, systems of coupled differential equations, and also implicit 
specification of the fc-derivative (as they are encountered in the fiow equation for a field dependent wave function). 

This paper is organized as follows: In Sec. |ll]we describe the numerical setup for FlowPy. in Sec. IIIII we describe 
how FlowPy is used, which parameters it takes and how they are specified. In I VII we discuss both simple and typical 
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examples of differential equations in the study of the Functional Renormalization Group. We give detailed examples 
for Python scripts that demonstrate how FlowPy can be used in practice. Our aim is to demonstrate how FlowPy is 
used, thus we just take flow equations from the literature without any derivation. As FlowPy is a general solver for 
partial differential equations we will denote the flow functions with f{k,x) most of the time, unless we consider some 
special physical system. 

II. A NUMERICAL FRAMEWORK FOR RG FLOWS 

The main design choices when building a framework for renormalization group flows are about how to support a 
reasonably large class of interesting problems while keeping code complexity manageable. 

One of the most challenging aspects of functional renormalization group problems is that the rate of change of the 
flow function f(k,x) can, for each value of the scale parameter k and at each point x, receive contributions from 
all other points x' . Hence, we are dealing with non-local partial integro-differential equations, or coupled systems of 
such equations, which may furthermore contain higher derivatives with respect to the coordinate x, and potentially 
be given in implicit form only. 

Conceptually, a numerical approach to problems of this nature involves numerical ODE solving (after discretization 
of the X range, as discussed below), numerical interpolation, numerical integration, as well as numerical differentiation. 
While an advanced numerical approach to such problems would perhaps develop a sophisticated combined (adaptive) 
discretization scheme that handles these different aspects in an unified fashion, we here contend ourselves with 
combining functionality from readily available libraries (specifically, functions from 0DEPACK[25'], FITPACK[2^, 
and QUAD PACK 27]) to solve the first three of the aforementioned tasks. Rather than working with these libraries 
directly, we use wrappers available in Scientific Python [28| . as the flexibility provided by the Python programming 
language is very helpful for addressing some subtle aspects of the task. 

The FlowPy package's objective is to numerically solve equations (resp. equation systems) of this type through 
discretization of the x-range with x denoting for example a field variable or a momentum. Choosing a number of 
support points Xj,j € {1, 2, . . . , jmax}, a PDF for /(fc, x) of the form 

^f{k,x)^G[k;x;f;d,f;dlf;...] (2) 

(where G is some functional) gets turned into a set of coupled ODEs of the form d/dkfj{k) = Gj[k; f; f; f"; . . .] 
with fj{k) = f{k,Xj) and /' being the vector dxf{k;xj) etc.. If the functional G involves integration over the 
position parameter (as it often does), the computational effort needed to evaluate the right hand side once (in order 
to numerically integrate the ODE) will grow quadratically with the number of support points. At the time of the first 
release of FlowPy, choosing the number of support points jmax somewhere between 20 and 100 and using a geometric 
distribution for the Xj seems an appropriate choice for many problems. 

An earlier prototype of FlowPy, which was used to do the calculations underlying (29| . was only concerned with 
performing the ODE integrations after x-discretization and required writing low level C code to specify the flow 
equation right hand side. It was soon found that this was a fairly tedious procedure that in particular needed some 
quite speciflc computing expertise beyond what may reasonably be expected from physicists wanting to solve RG 
flow equations. For this reason, the version of FlowPy described here was extended with an equation parser and 
code generator that automatically translates flow equation specifications to machine code and then loads this for fast 
execution. This 'equation compiler' is now the largest component of the FlowPy package . A similar program package 
to automate the calculations of Dyson-Schwinger equations was presented in [30|. 

III. USING FLOWPY 

FlowPy needs the following packages to be pre-installed: 

• Python {2.x with a; > 6) 

• Scientific Python (SciPy) 

• The 'NumPy' package for specialized numerical arrays 

• The 'python-simpleparse' parser generator package 
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• AC compiler such as gee (the default compiler used by FlowPy), as well as the Python library header files 
(especially Python. h). 

If the Python extension packages are installed into a location not normally searched by Python, the environ- 
ment variable PYTHDNPATH must be configured to include the installation-specific Python extension module path (e.g. 
$H0ME/lib/python2 . 7/site-packages). 

The source code zip file that accompanies this article can be downloaded either together with the arXiv preprint 
source from|http: / /arxiv.org/e-print / ARTICLE- ID, or through the 'ancilliary files' link on arXiv. This zip file contains 
the FlowPy code, plus html documentation, license and installation information, and some examples. 

Python header files are required by FlowPy as it will use the C compiler to translate user-specified equations to 
a a compiled Python extension module. This FlowPy-generated C module refers to a number of low-level Python 
definitions. On Linux systems, these header files usually come in a package named python2 . 7-dev or similar. 

If these packages are installed and Python is configured correctly, running python test.py in a console in the 
project folder performs a test. If it is successful the output is similar to: 



Users are strongly advised to perform this check in order to ensure FlowPy has been installed and set up correctly 
before using it. 

A complete example showing basic use of FlowPy is provided in the file demo . py. This shows how to compute a 
renormalization group fiow from k = 10^ to A: = 10"'^ with 20 intermediate fc-steps in geometric distribution for the 
flow equation 

-TT 

with 

C-{pi,p2,(t>) ^ pI+pI -2pip2COS(j> and C+{pi,p2, (f)) = pi + pi + 2pip2 cos (f> (4) 
The demo . py Python code is shown in Listing [TJ 

Listing 1: The demo.py Python code 



# demo . py 

from FlowPy import grange, flowproblem, make_f low_logger , make_lhs_iterator 

# Note: triple-quote """ bounds multi-line Python strings 
eqns= 

d/dk f(k,x) = 
integral[d phi from -pi to pi, 
d q from to k] 
(Cminus(x~2-k"2,q,phi)/Cplus(x~2-k~2,q,phi))* 
(1- d/dk f (k,x)) ; 

FLOWSTART f(k,x)=l; 

# === Helper functions === 

Cminus (pi , p2 , alpha) = pi ~2+p2 ~2-2*pl*p2*cos ( alpha) ; 
Cplus (pi , p2 , alpha) = pi ~ 2 + p2 ~ 2 + 2*pl *p2 * cos ( alpha ) ; 

# === Constants === 

pi = 3.241592653587983; 
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xs_plus=grange (0 . 1 , 100 , 10) # geometrically subdivided interval 
xs=[-xj for xj in xs_plus] + [0 . 0] +xs_plus # x-discretization 

f p=f lowproblem ( " demo_problem" , # problem/project name 

xs , # position discretization 

equations=eqns , # the equations 

# Perform 20 k-steps geometrically distributed 

# between 10~5 and 10~(-3): 
ks= grange (le5,le-3,20) , 

# Log flow data to file: 

log_state=make_flow_logger("demo_problem.flow"), 
# 

# For implicit flow equations that involve 

# the LHS expression (here d/dk Z2(k,p)) 

# on the RHS : how to handle iterative determination 

# of a self -consistent RHS value: 
# 

decide_iterate=make_lhs_iterator(loops=0) , 
# 

# Discretization Parameters below 

# 

# Step size for higher-order numerical differentiation: 
eps_dif f =le-4 , # this also is the default. 
diff_ord=4, # this also is the default. 

interpolation_kind=4 , # 4th order interp (also default). 
verbose=2, # verbose=0 turns off extra flow debug messages 
) 

fp . f lowO 

# After this, data will be in the file " demo_problem . f low" . 

The FlowPy package provides the f lowproblem class as well as some auxiliary functions such as grange (to produce a 
geometric distribution of numbers in a given range) and linrange (to produce an arithmetic distribution of numbers in 
a given range), make_f low_logger (to produce an object that can be used as log_state parameter to the f lowproblem 
constructor function) and make_lhs_iterator (to produce a decision function that can be used as decide_iterate 
parameter). The only relevant method of the flow problem class is the obj .flow() method that executes the solution 
of the problem. 

In order to solve the flow equation, a f lowproblem object is created. The default values for parameters are given 
in the following: 

FlowPy . f lowproblem (problem_name , 

xs , 

equations , 

ks= grange (le5,le-3,20) , 

log_state=None , 
eps_dif f =le-4 , 
dif f _ord=4 , 
Interpol at ion_kind=4 , 

decide_iterate=make_lhs_iterator(eps_abs=le-8), 

verbose=0 , 

) 



The mandatory parameters to FlowPy. flowproblemO are: 
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problem_n£mie The problem name - this is also used to create a directory that will contain 
machine-generated low level code which has been produced from the user spec- 
ified equations 

xs a list of support points for x-discretization 

equations the flow equations, as a (usually long) string. These will be handed to the 
parser and code generator. 

Other parameters that can be used for writing a logfile, writing to standard output in order to allow immediate 
supervision of the flow, interpolation order, steps for the k integration etc. are: 

ks /c-steps for the flow 

log_state function with call signature f (ff_ncmies,xs,k,ff_ys,ff_ydots) that 

logs a flow state. Default None means: Do nothing. ff_ys and 
ff_ydots are lists of arrays, one for each flow function, containing 
the values at the discretized support points xs for the given value of k. 
eps_dif f Step size for numerical differentiation in right hand side expressions 

dif f _ord Numerical differentiation will be done in a way that is correct to this 

given order 

interpolation_kind Specifies interpolation to be used (as a parameter to 
scipy . interpolate . interpld). 

decide_iterate A decision function f mapping f (k, history) to True / False, True 
meaning "do another iteration to determine LHS d/ dk" (cf . documen- 
tation of inake_lhs_iterator) 

verbose Determine verbosity level for reporting 

cc_call Template pattern for calling the C compiler. (Default: 

gcc 2>$F.cc_log -I /usr/include/python$V -fPIC 
-shared -o $F.so $F.c) 

The flow equations are defined as a string in between delimiters " " " and consist of a collection of definitions, each 
ending in a semicolon ' ; '. There are four possible kinds of definitions: 

• Constant definitions such as c=100; 

• Helper functions such as hb(vdiff ,k,z)=vdiff+k*z~2; 

• Starting conditions for the flow equation such as FLOWSTART f (k,x)=l;. 

• Flow equations such as d/dk f (k,x) = [rhs] ; The first argument in the functions always has to be the RG- 
parameter fc, the second can be a field variable, a momentum etc. 

There are two kinds of boundary conditions: For each fiow function, we have to specify the values at the starting 
value of k. This is done with a FLOWSTART definition. At the boundaries of the p-discretization range Smin/maxj the 
flow will, for every fiow function / and every value of fc, be clamped to (d/dk) f{k,x^in /ma,^) = 0. It is planned to 
drop this restriction in a future release of FlowPy to allow for more flexible boundary conditions. 

The order of definitions does not matter because FlowPy will sort them automatically into the right order, and 
hence they can be written down in the way that best describes the problem. However, there are of course some obvious 
constraints: Definitions of constants and auxiliary functions may involve auxiliary functions or other constants, but 
of course there must not be circular dependencies in these definitions. The compiler will detect and report circularity 
if this rule is violated. Also, while flow functions may be implicit i.e. the right hand side depends on the left hand 
side, and may involve auxiliary functions, fiow functions must not be used on the right hand side of the definition of 
an auxiliary function or constant. 

Flow equations can be coupled, given implicitly, can contain at most two integrations and can contain derivatives 
(also of higher order) with respect to the second argument. FlowPy will report an error if fiow equations contain 
derivatives too high for the interpolation scheme used. A number of basic examples are discussed in Section IVII 

IV. FLOWPY'S INTERNAL MECHANICS 

In order to use FlowPy to full effect, it is useful to have a basic mental model of its internal design. When a 
f lowproblem object is created, this will use FlowPy 's built-in compiler to translate the user-specified flow problem 
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equations to C code that subsequently gets compiled (calling the external C compiler) to a Python extension module. 
The C code, as well as the Python module object code and a compilation logfile will be placed in a special directory 
that is created with the name of the user-provided flow problem. These low-level files are named f lowpycoreN. *, 
where N starts at 1 and is increased by one for every new problem produced by the same Python process. (This 
is needed to inform the Python module import mechanism that different f lowpycoreN . so objects specify different 
problems.) When running multiple FlowPy processes simultaneously on a cluster from within a cluster-wide visible 
directory, it is strongly advised to ensure that the problem_name parameter contains a task-id so that different tasks 
running on different computer nodes use different auxiliary directories for their low-level modules. 

The grammar specifying flow problems is defined in the FlowPy file FlowPy^rammar . py; the parser that produces 
the parse tree from a problem specification is given in FlowPy_parser .py. The file FlowPy _builtins .py contains 
definitions of built-in FlowPy functions (such as exp, cos, etc.), as well as the skeleton for the C code of the Python 
extension module to be generated. The generated Python module will contain code that initializes all constants, defines 
user-specified auxiliary functions, parameter-dependent integration boundaries, initial conditions for fiow functions, 
as well as flow-function right-hand sides. When the evaluation of a flow function right-hand side needs to access an 
interpolated value of the flow function, or any of its derivatives, the C code will perform a callback into Python, 
evaluating a generic interpolating function that was provided to it by FlowPy. A considerable amount of code magic 
is hidden in these Python-defined interpolating functions that also can provide interpolated values for derivatives. 
Typically, these are complex closures involving various SciPy functions that interface FORTRAN code. 



V. FORMAT OF OUTPUT FILES 

In Listing [5] we show as an example part of a logfile from Sec. I VI B] provided by FlowPy. 



Listing 2: Sample output provided by FlowPy 



# === 


FLDWPYl 


. 


LOG 


FILE === 












# xs = [ 


-1 






> 


-0 


909091 


+ 




+ 


909091 


+ 1] 


# f lowf uns= [ ' 


E' , 


' lambda ' , ' omega 


'] 










# k ; 


nr. 


.flowfun ; 


ys/ydots=0/l ; 


y s / ydot s [0] ; 


ys/ydots [1] ; 




+ 100000 








+ 


+ 




+ 




+ 




+ 


+ 100000 


1 





+ 1 


+ 1 




+ 1 




+ 1 




+ 1 


+ 100000 


2 





+ 1 


+ 1 




+ 1 




+ 1 




+ 1 


+ 100000 





1 


+ 


-3 


1831e-ll 


-3 


1831e-ll 


-3 


1831e-ll 


+ 


+ 100000 


1 


1 


+ 


+ 1 


90986e-20 


+ 1 


90986e-20 


+ 1 


90986e-20 


+ 


+ 100000 


2 


1 


+ 


-3 


1831e-ll 


-3 


1831e-ll 


-3 


1831e-ll 


+ 


+ 39810 


. 7 








+ 


+ 4 


90063e-06 


+ 4 


90063e-06 


+ 4 


90063e-06 


+ 


+ 39810 


. 7 


1 





+ 1 


+ 1 




+ 1 




+ 1 




+ 1 


+ 39810 


. 7 


2 





+ 1 


+ 1 




+ 1 




+ 1 




+ 1 


+ 39810 


. 7 





1 


+ 


-2 


00841e-10 


-2 


00841e-10 


-2 


00841e-10 


+ 


+ 39810 


. 7 


1 


1 


+ 


+ 7 


60329e-19 


+ 7 


60329e-19 


+ 7 


60329e-19 


+ 


+ 39810 


. 7 


2 


1 


+ 


-2 


0084e-10 


-2 


0084e-10 


-2 


0084e-10 


+ 


+ 15848 


. 9 








+ 


+ 1 


70471e-05 


+ 1 


70471e-05 


+ 1 


70471e-05 


+ 


+ 15848 


. 9 


1 





+ 1 


+ 1 




+ 1 




+ 1 




+ 1 


+ 15848 


. 9 


2 





+ 1 


+ 1 


00002 


+ 1 


00002 


+ 1 


00002 


+ 1 


+ 15848 


. 9 





1 


+ 


-1 


26724e-09 


-1 


267246-09 


-1 


26724e-09 


+ 


+ 15848 


. 9 


1 


1 


+ 


+ 3 


02692e-17 


+ 3 


02692e-17 


+ 3 


02692e-17 


+ 


+ 15848 


. 9 


2 


1 


+ 


-1 


26721e-09 


-1 


26721e-09 


-1 


26721e-09 


+ 


+6309 . 


57 








+ 


+ 4 


74848e-05 


+ 4 


74848e-05 


+ 4 


74848e-05 


+ 


+6309 . 


57 


1 





+ 1 


+ 1 




+ 1 




+ 1 




+ 1 


+6309 . 


57 


2 





+ 1 


+ 1 


00005 


+ 1 


00005 


+ 1 


00005 


+ 1 


+6309 . 


57 





1 


+ 


-7 


99596e-09 


-7 


99596e-09 


-7 


99596e-09 


+ 


+6309 . 


57 


1 


1 


+ 


+ 1 


20504e-15 


+ 1 


20504e-15 


+ 1 


20504e-15 


+ 


+6309 . 


57 


2 


1 


+ 


-7 


99558e-09 


-7 


99558e-09 


-7 


99558e-09 


+ 



On the second line, starting with # xs are the discretization points, the third line gives the functions for which the 
fiow equations were solved. Starting from line five, the results from the solution of the fiow equation are listed. The 
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first column denotes the k values, the second column is a number corresponding to the respective function, E, lambda 
and omega in the logfile above. and 1 in the third column stand for the function and its derivative respectively. 
The following columns display the values of the functions evaluated at the discretization points. As can be seen, the 
values at the boundary are fixed throughout the fiow. 

VI. EXAMPLES 

In this section we will show how typical examples for flow equations are solved with FlowPy and the results are 
compared with solutions obtained from Mathematica and SciPy, in order to establish that FlowPy solves these 
problems correctly. In the subsequent sections we will also display examples of Python scripts which demonstrate the 
specification of the flow equations and the parameters that FlowPy takes. 

A. Simple examples 

The first five examples are devoted to the solution of simple differential equations, most of which have an analytic 
solution. The purpose of these first examples, which are contained in the test.py test file, is to demonstrate in a 
readily verifiable way that FlowPy correctly works as claimed. 

1. Constant growth 
As a first example we solve the differential equation 

^/(fc:^) = -l (5) 

describing a constant growth with the fiow parameter k lying between k £ [A, ko] with A = 110. The function f{k, x) 
is specified at the scale A and the boundary conditions for Xstart = and a;end = 10 are chosen in the following way 
at fc = A and x — a;start /end ■ 

/(A, x) = X, f{k, Xstart) = /(A, Xstart) .f{k, Xcnd) = /(A, Xond) (6) 

The analytic solution of this equation with the above starting condition is 

f{k,x) ^ A~ k + x ^ no- k + X. (7) 

The corresponding FlowPy problem specification is (note that d/dkf is specified as being negative as we are fiowing 
from large to small k): 

f lowpr oblem ( " example 1 " , 

xs=[float(n) for n in xrange(ll)] , 

equat i ons = " " " 

d/dk f(k,c)=-l; 

FLDWSTART f(k,x)=x; 
II II II 

log_state=logger 1 , 
ks= grange (110,10,1) , 

decide_iterate=make_lhs_iterator (loops=0) ) 

Solving this fiow equation with FlowPy yields the correct solution as can be easily checked with eq. [T] 



X 





123456789 


10 


/(fc = 10,x) 





101 102 103 104 105 106 107 108 109 


10 



As expected the function at the boundary is fixed to the values at k = A throughout the whole flow. 
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2. Exponential growth 

The second example is a differential equation describing exponential growth, 

^/(fc,x) = -0.01-/(fc,a:), (8) 
with A — 110, Xstart = and a;cnd — 10 and the values at the boundary fixed to the value at fc = A for all values of fc, 

/(A,a;)=X, /(fcjXstart) = /(A, a;start), /(fcjXond) = /(A,Xcnd). (9) 

The solution of this equation with the above starting condition is 

= a;e-°-°i('=-^) (10) 
In FlowPy the equations are specified in the following way: 



f=flowproblem(" example 2 " , 




xs=xs , 




equat ions = " " " 




d/dk f(k,x) = -0.01*f(k, 


x) ; 


FLOWSTART f(k,x) = x; 




It It It 




log_state=logger2 , 




ks = grange (110 , 10 , 1) , 




decide_iterate=make_lhs_ 


iter at or (loops=0) ) 



The numerical solution of this equation with FlowPy yields 



X 


1 2 


3 4 5 6 7 8 9 


10 


/(fc = 10,x) 


2.71828 5.43656 I 


5.15485 10.8731 13.5914 16.3097 19.028 21.7463 24.4645 


10 



which is in accordance with the analytic solution from eq. (|10p . Again, the values at the boundary do not flow. 

3. Implicit d/dkf{k, x) 

Especially when considering a problem additionally containing a wave function renormalization, the flow equation 
for the wave function renormalization is given implicitly. In order to demonstrate that FlowPy can also solve implicitly 
given functions, we discuss the differential equation given by 

^/(fc,x) = -1.0 + 0.5^/(fc,x) (11) 

with the starting condition 

fiA,x) = x (12) 
and a:start, a^end and A the same as above. The analytic solution to this equation is 

f{k,x) = 2A- 2k + x (13) 

The FlowPy code reads: 



f=flowproblem(" examples " , 
xs=xs , 

equat ions = " " " 

d/dk f(k,x) = -1 . 0+0 . 5*d/dk f(k,x); 

FLOWSTART f(k,x) = x; 
II II II 

log_state=logger3 , 
ks = grange (110 , 10 , 1) , 

decide_iterate=make_lhs_iterator (eps_abs=le-6) ) 



and solving numerically with FlowPy yields 

X 0123456789 10 



/(A; = 10,a;) 



201 202 203 204 205 206 207 208 209 10 
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4- Differential equation with integral on the right hand side 

Often in the study of the renormalization group equations, the right hand side is given as an integral in the 

momentum. FlowPy can handle also such a right hand side as is dc^inonstratcd in this section. The flow rate is 
constant for each point, but for each a;-value, we express the flow rate as an integral: 

^f{k,x) = - J\yy^ (14) 

The starting condition is chosen to be 

f{k = A,x)=x (15) 

and the solution to this equation is 

f{k,x) = ^{A-k)x^+x. (16) 

The 'flowproblem' takes the form 

fp=flowproblem(" example4" , 
xs=xs , 

equat ions = " " " 

d/dk f(k,x) = integral [dq from to x] -q~2; 

FLOWSTART f(k,x) = x; 
It It It 

log_state=logger4, 
ks = grange (110 , 10 , 1) , 

decide_iterate = inake_lhs_iterator(loops=0)) 
and the numerical solution with FlowPy is 



.)■ 


1 2 3 4 5 () 7 8 9 10 


/(fc = 10,x) 


34.3333 268.667 903 2137.33 4171.67 7206 11440.3 17074.7 24309 10 



Again, the boundary values are flxed throughout the flow. 

5. Heat equation 
As a more complicated example we solve a heat flow problem 

■^f{k,x) = -0.01 ■ f"{k,x), /(A; = A,x) = exp(-0.5(a;-5.0)V(2.02)) (17) 



with FlowPy. The flowproblem reads: 



fp = flowproblein(" examples " , 




xs = xs , 




equations=" " " 




d/dk f (k ,p) = -0 . 01*f ' ' 


(k,p) ; 


FLOWSTART f(k,p) = exp ( 


-0.5*(p-5.0)-2/(2.0-2)) ; 


It It It 




ks=ks , 




log_state=logger5 , 




decide_iterate=make_lhs 


_iterator(loops=0)) 



Values obtained with the odeint routine from SciPy and a somewhat less sophisticated discretization of the second 
derivative are given the table below. 
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Xs 


solution from SciPy 


solution from FlowPy 


deviation in % 





0.04394 


0.04394 


0.000 


0.1 


0.05730 


0.05731 


0.006 


0.2 


0.07068 


0.07069 


0.010 


0.3 


0.08410 


0.08411 


0.012 


0.4 


0.09756 


0.09757 


0.013 


0.5 


0.11109 


0.11110 


0.014 


0.6 


0.12468 


0.12470 


0.015 


0.7 


0.13837 


0.13839 


0.014 


0.8 


0.15215 


0.15217 


0.013 


0.9 


0.16603 


0.16605 


0.013 


1 


0.18001 


0.18003 


0.012 


1.1 


0.19410 


0.19412 


0.011 


1.2 


0.20829 


0.20831 


0.010 


1.3 


0.22258 


0.22260 


0.009 


1.4 


0.23697 


0.23699 


0.008 


1.5 


0.25144 


0.25146 


0.006 


1.6 


0.26599 


0.26600 


0.005 


1.7 


0.2806 


0.28061 


0.004 


1.8 


0.29526 


0.29526 


0.002 


1.9 


0.30995 


0.30995 


0.001 


2 


0.32465 


0.32465 


0.000 



If the number of discretization points for x is large enough, the deviation will become arbitrarily small. (Due to the 
simplistic discretization method used by the SciPy ODE solver based code, the FlowPy results are likely closer to the 
exact solution here.) 

Now that we have convinced ourselves that FlowPy solves these trivial examples correctly, we turn to typical 
equations that are encountered in the study of functional rcnormalization group equations. 



B. System of coupled ordinary differential equations 

In this section we solve the flow equations for the anharmonic oscillator. This will give us an example for a system 
of coupled flow equations. The flow equation for the effective potential, using the regulator Rk = (fc^ —p'^)9{k'^ —p'^), 
is given by [(3] 

where prime denotes the derivative with respect to x. In order to obtain a system of ordinary differential equations 
we make a polynomial approximation for Vk{x) — Ek + WfcX^ + A^a;'^. This yields 

d ~ _l / fc2 \ d _ 1 fe^ d _ 6 ^2 

'dk^''-n\WT^k'^)' dk"^' ' 'nWT^^'' dk^' - nWT^^"- 

Because FlowPy by design assumes flow equations to be partial differential equations involving one extra parameter 
beyond the rcnormalization scale fc, we have to use a small workaround if we want to solve a system of ordinary 
differential equations. We have to insert an artificial x dependence in all couplings and solve this system of pseudo 
partial differential equations. The values for couplings is constant for all x G [xstart+i, 2;cnd-i]- 

A python script to solve these differential equations is shown in Listing [3l 

Listing 3: An example for a system of ordinary differential equations 



# anharmonic_oscillator . py 
import FlowPy 
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# specifying the flow equations: 
aho_eqns= " " " 

d/dk E(k,x)=(l/pi)*(k~2/(k"2+omega(k,x))-l) ; 

d/dk omega (k , x ) = -k ~ 2* lambda (k,x)/(pi*(k~2 + omega(k,x))~2) ; 

d/dk lambda (k,x)=(6/pi) *k~2/ (k~2+omega (k , x) ) ' 3* lambda (k , x ) ~2 ; 

♦specifying the starting conditions: 
FLOWSTART E(k,x)=0; 

FLOWSTART omega (k , x )= 1 ; 
FLOWSTART lambda (k , x) =10 ; 

# === Parameters === 
pi=3. 141592653589793; 



# specifying the values for x: 
xs_plus= [1 . l**n for n in range(3)] 

xs=[-x for x in reversed ( xs_plus )]+ [0 . 0] +xs_plus 

logger = FlowPy.make_flow_logger(filename = " anharmoni c_oscillat or . f low " ) 

# solving the flow equation: 

fp = FlowPy . flowproblemC anharmoni c_oscillat or " , 

xs=xs , 

equations=aho_eqns , 
log_state=logger , 

decide_iterate=FlowPy.make_lhs_iterator(loops=0), 

verbose=3 

) 

fp . f lowO 

The results for Ek obtained with FlowPy for different values of Aa are displayed in the table below and coincide 
to the third digit with those obtained for example with the DSolve routine from Mathematica 5. This deviation is 
mostly due to the different discretizations in k. 



Aa 


FlowPy 


Mathematica 





0.49968 


0.49998 


0.1 


0.50276 


0.50306 


0.2 


0.505756 


0.506075 


0.3 


0.508676 


0.508994 


0.4 


0.511524 


0.511842 


0.5 


0.514306 


0.514624 


1 


0.527365 


0.527683 


1.5 


0.53928 


0.539598 


2 


0.550305 


0.550623 


2.5 


0.56061 


0.560928 


3 


0.570315 


0.570634 


4 


0.588268 


0.588586 


5 


0.604666 


0.604985 


6 


0.619835 


0.620154 


7 


0.633999 


0.634317 


8 


0.647321 


0.647639 


9 


0.659924 


0.660242 


10 


0.671905 


0.672223 
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C. Field dependent flow equations 

As an example of a field dependent flow equation we consider supersymmetric quantum mechanics. This model 
describes an anharmonic oscillator coupled in a supersymmetric way to fermions. The flow equation is derived for the 
superpotential W{(j)) which enters in the scalar potential a.s V = ^W'{(j)). The flow equation for the superpotential 
reads [15] : 

A python script to solve this equation is given in Listing 2] 

Listing 4: An example for a field dependent flow equation 



# SUSY_qM.py 
import FlowPy 




SUSY_qM_eqns= 

d/dk V(k,phi)=-V' ' (k , phi ) /4/ ( V ' (k , phi ) +k ) 
FLOWSTART V(k,phi) = l+m*phi+g*phi ~2+a*phi ~3 ; 


-2; 


# === Parameters === 
m = l ; 
g = 0; 
a = l ; 

It M ir 




#change n_start and n_end to change the phi range 

n_start=-2 

n_end=2 

n_steps=10 #increase the number of sampling points to improve the resolution 


xs=FlowPy . linrange(n_start , n_end , n_steps) 




logger=FlowPy.make_flow_logger(filename="SUSY_QM. 
fp=FlowPy .flowproblem("flowpy_SUSY_QM", 

xs = xs , 

equations = SUSY_QM_eqns , 
log_state = logger , 

decide_iterate=FlowPy . make_lhs_ 

verbose=2 , 

) 


flow") 

it er at or ( loops =0) , 


f p . f low ( ) 





The solution to the flow equation is shown in Fig.[T] Note that this picture was done with a resolution of n_step=40. 
The results are consistent with those obtained with the NDSolve routine with Mathematica 7 in llSl] 



D. System of two coupled field dependent flow equations 

As an example for a system of two coupled flow equations we consider supersymmetric quantum mechanics with 
field dependent wave function renormalization. 

As in the previous section, the flow equations with a field dependent wave function renormalization are discussed 
in [ll]. They read 
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FIG. 1: The scalar potential V = ^W'{4>) before and after the flow 



where we have introduced the abbreviations 



AA= (l + /c9fc)Zfe(0)2 and V ^ W" {(j)) + kZk{(p)^ 



(22) 



The python script to solve these equations is given in Listing [S) The solutions are shown in Fig [2] This picture was 
created with ii_steps=40. 

Listing 5: An example for a system of two coupled differential equations 



# SUSY_QM_wavef unction . py 
import FlowPy 

d/dk V(k,phi)=-V''(k,phi)*(Z(k, phi )~2+2*k*Z(k,phi)*d/dk Z(k,phi)) 

/4/ (V '(k,phi)+k*Z(k,phi)~2)~2; 
d/dk Z(k,phi)=l/Z(k,phi)*(Z(k, phi ) '2+2*k*Z (k , phi) * d/dk Z (k , phi) ) 

/4/(V'(k,phi)+k*Z(k,phi)'2)~2 

*(4*Z'(k,phi)*V''(k,phi)/(V'(k,phi)+k*Z(k,phi)~2) 
-Z' ' (k,phi)*Z(k,phi)-Z' (k , phi ) *Z ' (k , phi ) 

-3*Z (k , phi) '2*V ''(k,phi)'2/4/(V'(k,phi)+k*Z(k,phi)"2)"2); 

FLOWSTART V(k,phi) = l+phi+a*phi '2+b*phi ~3+c*phi ~4 ; 
FLOWSTART Z(k,phi) = 1; 



#==== Parameters ==== 
a = l ; 
b = l ; 
c = l : 



#change n_start and n_end to change the phi range 

n_start=-2 

n_end=2 

n_steps=10 #increase the number of sampling points to improve the resolution 
xs=FlowPy . linrange(n_start , n_end , n_steps) 

logger = FlowPy. make_flow_logger(filename = " SUSY_QM_wavef unct i on . flow") 
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FIG. 2: left panel: Scalar potential V = , right panel: field dependent wave function 



fp=FlowPy . flowproblem("flowpy_SUSY_QM_wavef unction" , 

xs=xs , 

equations=SuSy_eqns , 
log_state=logger , 

decide_iterate=FlowPy . make_lhs_iterator (loops=0) , 

verbose=2 , 

) 

fp . f low () 



E. Momentum-dependent flow equations 

FlowPy can also solve flow equations with full momentum dependence. In this example we calculate the momentum 
dependent wave- function renormalization for the supersymmetric N = 2 Wess-Zumino Model in two dimensions. This 
model is thoroughly discussed in [2^ where the flow equation is derived. The flow equations read 

d,ZUp) = -I6g^ J ^J^J^^^Zl (,) ZU\P~<1\) du [kZl (,)) , (23) 

where we have introduced the abbreviation 

N{q) = {q^Zt (q) + {kZl (q) + mf) . (24) 

q and p denote two dimensional vectors, whereas q = |q| and p — \p\ respectively. 

A python script to solve this equation with m = 1 and g = 0.3 is given in Listing |6l Note that due to the complex 
structure of the flow equations this script takes a long time to calculate the wave function even with a very low 
resolution of five sampling points. The result is shown in Fig. [S] The calculation was done with n_steps=30. 

Listing 6: An example for a momentum dependent wave function 

# moment um_dependent_wave function . py 
import FlowPy 
import math 

wavef unction_eqns=" " " 
d/dk Z(k,p) = 
integral[d phi from -pi to pi, 
d q from to k] 

(-16*g"2)/4/pi"2*(k*Z(k,q)+m)=f=Z(k,q)*Z(k,Cminus(p,q,phi)) 
*(l+k*d/dk Z(k,q))/(q*q*Z(k,q)*Z(k,q)+(k*Z(k,q)+m)*(k*Z(k,q)+m)) 
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/ ( Cminus (p,q,phi)* Cminus (p,q,phi)*Z(k, Cminus (p,q,phi)) 
*Z (k , Cminus (p,q,phi))+(k*Z(k, Cminus (p ,q,phi) )+m) 
* (k*Z (k , Cminus (p , q , phi ) ) +m) ) 

/(q*q*Z(k,q)*Z(k,q)+(k*Z(k,q)+m)*(k*Z(k,q)+m))-2; 
FLOWSTART Z(k,p) =1; 

# === Helper functions === 

Cminus (pi , p2 , alpha) = pi '2+p2 ~2-2*pl*p2*cos ( alpha) ; 

# === Parameters === 
m = 1.0; 

g=0.3; 

pi=3. 141592653589793; 
II II II 

#change n_start and n_end to change the p_range 

n_start=-20 

n_end=10 

n_steps=5 ttincrease this number of sampling points to improve the resolution 

# logarithmic -equidistant distribution of the momentum: 
qs_log=FlowPy. linrange(n_start , n_end , n_steps) 

qs= [-math . exp (q) for q in qs_log] + [0 . 0] + [math . exp (q) for q in qs_log] 

# specifing the logfile: 

logger = FlowPy . make_flow_logger(filename = " moment um_dep endent _wavefunct ion . flow") 

# solving the flow equation: 

fp=FlowPy . flowproblem("flowpy_wavefunction" , 

xs=qs , 

equations=wavef unction_eqns , 
log_state=logger , 

decide_iterate=FlowPy . make_lhs_iterator (loops=0) , 
verbose=2) 

f p . f low ( ) 



F. Parametric Example 

In some situations, one would want to perform a number of renormalization group flow calculations that differ only 
in the choice of some parameters. While this can be achieved in a reasonably straightforward way with FlowPy via 
Python scripting, the FlowPy package contains a few extra definitions to make this process more convenient for the 
user, which are explained by means of an example in Listing [71 

Listing 7: Example for a flow equation with two parameters 

# parametric . py 

from FlowPy import flowproblem, flowparams, make_lhs_iterator , \ 
grange, make_f low_logger 

runs_todo= [flowparams (mu=val_mu ,T=val_T) 

for val_mu in [1.0,1.5,2.0,2.5] 
for val_T in [0.1,0.2,0.5,1.0]] 

xs=[float(n) for n in range(ll)] 
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1.004 




-10 -5 5 10 

ln(p/m) 

FIG. 3: Momentum dependent wave function renormalization 



for params in runs_todo: 
f = f lowpr oblem ( " mu_T " , 
xs=xs , 

equat i ons = " " " 
d/dk F(k,p) = -T; 
FLDWSTART F(k,p) = mu ; 
" " "-Hparams . def s () , 
ks = grange (110 , 10 ,5) , 

decide_iterate = make_lhs_iterator (loops=0) , 
log_state=make_flow_logger (\ 

params . subs _short ( " FLOW__mu=${mu}_T = ${T>__ . f low " ) ) , 
verbose=2) 

f . f lowO 



The key concept here is the f lowparsuns "FlowPy parameters" object. This represents a collection of parameter 
choices that can easily be mapped to strings in various contexts, e. g. to generate systematic output filenames. 
A f lowparamis object is created by specifying explicitly all "parameter name = numerical value" associations in a 
function call of the form 

params = f lowparams (mu=l . 5 , alpha=3.0, T=0.5) 

A potentially relevant limitation is that special keywords in the Python programming language (such as e.g. def, 
raise, lambda, else, etc.) can not be chosen as parameter names. (There is a way to avoid this problem in Python 
by using a different function call form, should this really turn out to be a problem.) 

If paramis is a FlowPy parameter object, e.g. defined as in the example in Fig. [7], then params. def s() will produce 
a string containing valid FlowPy constant definitions that introduce these parameters. One would typically want to 
append this to the (parameter-dependent) flow equation definitions by using Python "string addition" as shown in the 
example. Another useful snippet is parcmis.subs_sh.ort (pattern) which will regard pattern as a template string 
on which parameter substitution has to be performed according to the rules of Python's string . Template class, i.e. 
${xyz} will be substituted by the pretty-printed numerical value of the parameter xyz. (Detailed rules can be found 
in the Python documentation of the string. Template mechanism.) One may want to use this to automatically 
generate filenames, as in the example shown. 

A very useful Python feature one might want to employ in such situations is special syntax available to define lists 
from cartesian products. In Python, one can e.g. write 

[(x,y) for X in ranged, 5) for y in ranged, x)] 
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to obtain the list 

[(2, 1), (3, 1), (3, 2), (4, 1), (4, 2), (4, 3)] 

This syntax, which is an adoption of a very similar feature that probably was first made popular through the Haskell 
programming language, can be used to concisely specify fairly complex constructions, such as 

from fractions import gcd 

[(x,y,x*y) for x in range (5) for y in range (5) if gcd (x , y ) ==1] 
which gives the list 



[(0, 


1, 


0) , 


(1, 


0, 


0) , 


(1, 


1, 


1) , 


(1, 


2, 


2) , 


(1, 


3, 


3) , 


(1, 


4, 


4) , 


(2, 


1, 


2) , 


(2, 


3, 


6) , 


(3, 


1, 


3) , 


(3, 


2, 


6) , 


(3, 


4, 


12) , 


(4, 


1, 


4) , 


(4, 


3, 


12)] 



Users of FlowPy might find the basic form of this construct useful to define cartesian products of parameter choices, 
as in the example provided in this section. 

Evidently, if the approach shown here is used to perform N different numerical renormalization group flow cal- 
culations, then TV different machine code files will be produced automatically by FlowPy and loaded into Python. 
Although this process could be improved conceptually, this would presumably be only worthwhile for TV ^ 100. This 
is unlikely to be relevant for typical applications, as problems of typical complexity are expected to be partitioned into 
collections of < 100 before being submitted to a computing cluster anyway. When performing multiple FlowPy 
calculations on a cluster, one should ensure that each cluster job chooses a different name when defining a flowprob- 
lem. Otherwise, accidental unintended sharing of the directory used to produce machine code may result, and some 
calculations may fail or, even worse, end up using a wrong set of parameters and equations. Typically, an unique 
job ID will be taken from the program's argument list, from environment variables, or generated semi-randomly using 
the time and process id. Python makes this information available via the sys.argv variable containing program 
arguments (the sys module has to be imported first), the os.getenvO function that behaves like the corresponding 
C library function, the time. time () and related functions, and the random module. 

VII. CURRENT LIMITATIONS 

In its present form, the FlowPy package has a number of limitations. None of these actually are fundamental - they 
all could be overcome with some dedicated effort on the side of the FlowPy authors. Concerning the question how 
relevant these are, and where to focus further development effort on, the authors seek input from the renormalization 
group flow community. 

• FlowPy assumes every flow function to depend on precisely one extra parameter beyond the scale parameter k. 

• The right-hand side of flow functions may involve zero to two integrals, and the integration boundaries may 
only depend on left-hand side parameters, i.e. the inner integration boundary cannot depend on the value of 
the outer integration parameter. Also, term notation presently requires integral specifications to directly follow 
the equals sign, i.e. all extra factors must be pulled under the integral. Expressions such as C J^^_^ da f{a) 
(with C a numerical constant) hence must be re-written as daCf{a). 

• Currently boundary conditions in the direction are fixed to the classical values of the effective potential. 
However, more flexible boundary conditions that can be specified by the user might be more appropriate for 
some problems. 

• The range of special functions made available to the user so far is fairly limited. 
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• Error reporting could be improved to make it easier to find typos in the specification of the flow equations. 

• FlowPy suffers from some minor bugs in the default parser definitions inherited from python-simpleparse. In 
particular, a number such as '5e-3' is not a valid constant - this must be given as '5.0e-3' instead. 

• This version of FlowPy does not support parallelization. However, re-introducing MPI support is planned. 

• There are many opportunities to make FlowPy more efficient by improving its internal design. This, however, 
will require substantial low level changes. 



VIII. OUTLOOK 

In this paper we present the numerical toolbox FlowPy which is able to solve many typical partial differential 
equations encountered studies of the functional renormalization group. We hope that it will prove to be useful in 
the application of functional renormalization group techniques and that it will facilitate these studies in providing a 
powerful numerical tool to solve the differential equations. 

We plan to include MPI support in the next version of FlowPy. Additionally it is planned to make the boundary 
conditions in the field variable more flexible such that it will be possible to use e.g. the first loop approximation. Also 
there is a lot of room for improvement in the performance of the parser. Here we present the first version of FlowPy 
and we have planned to release updates in the future. The intention of this paper is to make FlowPy available to a 
wider audience. Thus we are grateful for suggestions how to improve FlowPy. 
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